<<<<<<< HEAD ======= >>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85 <<<<<<< HEAD <<<<<<< HEAD

Modeling Relative Survival/Excess Mortality using stpm2 (Part 2):

Flexible parametric modelling using restricted cubic splines

Miguel Angel Luque-Fernandez

Statistical methods for population-based cancer survival analysis (CIBERSP & VICA), Spain

Note: Adapted with permission from Paul Dickman and Paul Lambert short course in Trevisso, Italy

=======

Modeling Relative Survival/Excess Mortality using stpm2 (Part 2):

Flexible parametric modelling using restricted cubic splines

Miguel Angel Luque-Fernandez

Statistical methods for population-based cancer survival analysis (CIBERSP & VICA), Spain

Note: Adapted with permission from Paul Dickman and Paul Lambert short course in Trevisso, Italy

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD

1. Calling Stata software from Jupyter Notebook

=======

1. Calling Stata software from Jupyter Notebook

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

import os
<<<<<<< HEAD
os.chdir('/Applications/stata/utilities/')
from pystata import config
config.init('mp')
=======
os.chdir('/Applications/stata/utilities/')
from pystata import config
config.init('mp')
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD

  ___  ____  ____  ____  ____ ®
=======
  ___  ____  ____  ____  ____ ®
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
 /__    /   ____/   /   ____/      17.0
___/   /   /___/   /   /___/       MP—Parallel Edition

 Statistics and Data Science       Copyright 1985-2021 StataCorp LLC
                                   StataCorp
                                   4905 Lakeway Drive
                                   College Station, Texas 77845 USA
                                   800-STATA-PC        https://www.stata.com
                                   979-696-4600        stata@stata.com

Stata license: Single-user 4-core  perpetual
Serial number: 501706307451
  Licensed to: MALF
               IBS.GRANADA

Notes:
      1. Unicode is supported; see help unicode_advice.
      2. More than 2 billion observations are allowed; see help obs_advice.
      3. Maximum number of variables is set to 5,000; see help set_maxvar.

<<<<<<< HEAD

2. Load Data, merge expected mortality and fit inital model

=======

2. Load Data, merge expected mortality and fit inital model

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

%%stata
cd /Users/MALF/Dropbox/CANSURV-slides-course/Granada/
use melanoma
stset surv_mm, failure(status=1 2) scale(12) exit(time 120.5) id(id)
gen _age = min(int(age + _t),99)
gen _year = int(yydx + _t)
sort _year sex _age
merge m:1 _year sex _age using popmort,  keep(match master)
tab agegrp, gen(agegrp)
gen female = sex == 2

/* Fit initial model */

stpm2, df(3) scale(hazard) bhazard(rate) 
predict h1, hazard per (1000) ci
predict s1, survival ci

<<<<<<< HEAD

. cd /Users/MALF/Dropbox/CANSURV-slides-course/Granada/
=======
. cd /Users/MALF/Dropbox/CANSURV-slides-course/Granada/
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
/Users/MALF/Dropbox/CANSURV-slides-course/Granada

. use melanoma
(Skin melanoma, diagnosed 1975-94, follow-up to 1995)

. stset surv_mm, failure(status=1 2) scale(12) exit(time 120.5) id(id)

Survival-time data settings

           ID variable: id
         Failure event: status==1 2
Observed time interval: (surv_mm[_n-1], surv_mm]
     Exit on or before: time 120.5
     Time for analysis: time/12

--------------------------------------------------------------------------
      7,775  total observations
          0  exclusions
--------------------------------------------------------------------------
      7,775  observations remaining, representing
      7,775  subjects
      2,777  failures in single-failure-per-subject data
 43,384.625  total analysis time at risk and under observation
                                                At risk from t =         0
                                     Earliest observed entry t =         0
                                          Last observed exit t =  10.04167

. gen _age = min(int(age + _t),99)

. gen _year = int(yydx + _t)

. sort _year sex _age

. merge m:1 _year sex _age using popmort,  keep(match master)

    Result                      Number of obs
    -----------------------------------------
    Not matched                             0
    Matched                             7,775  (_merge==3)
    -----------------------------------------

. tab agegrp, gen(agegrp)

   Age in 4 |
 categories |      Freq.     Percent        Cum.
------------+-----------------------------------
       0-44 |      2,046       26.32       26.32
      45-59 |      2,238       28.78       55.10
      60-74 |      2,280       29.32       84.42
        75+ |      1,211       15.58      100.00
------------+-----------------------------------
      Total |      7,775      100.00

. gen female = sex == 2

. 
. /* Fit initial model */
. 
. stpm2, df(3) scale(hazard) bhazard(rate) 

Iteration 0:   log likelihood = -8731.8053  
Iteration 1:   log likelihood = -8590.4493  
Iteration 2:   log likelihood = -8590.0259  
Iteration 3:   log likelihood = -8590.0249  
Iteration 4:   log likelihood = -8590.0249  

Log likelihood = -8590.0249                              Number of obs = 7,775

------------------------------------------------------------------------------
             | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
xb           |
       _rcs1 |   .8252308   .0249859    33.03   0.000     .7762595    .8742022
       _rcs2 |   .2110309   .0235833     8.95   0.000     .1648085    .2572534
       _rcs3 |   .0631928   .0109672     5.76   0.000     .0416974    .0846881
       _cons |  -1.813097   .0314253   -57.70   0.000     -1.87469   -1.751505
------------------------------------------------------------------------------

. predict h1, hazard per (1000) ci

. predict s1, survival ci

. 

<<<<<<< HEAD

3. Plot the predicted hazard and survival functions

=======

3. Plot the predicted hazard and survival functions

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

%%stata

twoway	(rarea h1_lci h1_uci _t, sort pstyle(ci)) ///
		(line h1 _t, sort), name(h1,replace) legend(off) ///
<<<<<<< HEAD
		xtitle("time since diagnosis") ///
		ytitle("excess mortality rate (per 1000 py's)")
						
twoway	(rarea s1_lci s1_uci _t, sort pstyle(ci)) ///
			(line s1 _t, sort), name(s1,replace) legend(off) ///
		xtitle("time since diagnosis") ///
		ytitle("relative survival")
=======
		xtitle("time since diagnosis") ///
		ytitle("excess mortality rate (per 1000 py's)")
						
twoway	(rarea s1_lci s1_uci _t, sort pstyle(ci)) ///
			(line s1 _t, sort), name(s1,replace) legend(off) ///
		xtitle("time since diagnosis") ///
		ytitle("relative survival")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD

. 
. twoway  (rarea h1_lci h1_uci _t, sort pstyle(ci)) ///
>                 (line h1 _t, sort), name(h1,replace) legend(off) ///
>                 xtitle("time since diagnosis") ///
>                 ytitle("excess mortality rate (per 1000 py's)")
=======
. 
. twoway  (rarea h1_lci h1_uci _t, sort pstyle(ci)) ///
>                 (line h1 _t, sort), name(h1,replace) legend(off) ///
>                 xtitle("time since diagnosis") ///
>                 ytitle("excess mortality rate (per 1000 py's)")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

.                                                 
. twoway  (rarea s1_lci s1_uci _t, sort pstyle(ci)) ///
>                         (line s1 _t, sort), name(s1,replace) legend(off) ///
<<<<<<< HEAD
>                 xtitle("time since diagnosis") ///
>                 ytitle("relative survival")
=======
>                 xtitle("time since diagnosis") ///
>                 ytitle("relative survival")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

. 

<<<<<<< HEAD

<<<<<<< HEAD

<<<<<<< HEAD

4. Compare fitted values for various df

=======

4. Compare fitted values for various df

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

%%stata

foreach df in 2 4 6 {
<<<<<<< HEAD
	stpm2, df(`df') scale(hazard) bhazard(rate)
	predict h_df`df', hazard 
	replace h_df`df' = h_df`df' * 1000
	predict s_df`df', survival 
	estimates store df`df'
=======
	stpm2, df(`df') scale(hazard) bhazard(rate)
	predict h_df`df', hazard 
	replace h_df`df' = h_df`df' * 1000
	predict s_df`df', survival 
	estimates store df`df'
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
}

/*Plot the excess hazard functions*/
twoway (line h_df2 h_df4 h_df6 _t, sort lcolor(red blue black))

/*Plot the survival functions*/
twoway (line s_df2 s_df4 s_df6 _t, sort lcolor(red blue black)) 


estimates stats df2 df4 df6, n(2773)

<<<<<<< HEAD

. 
. foreach df in 2 4 6 {
  2.         stpm2, df(`df') scale(hazard) bhazard(rate)
  3.         predict h_df`df', hazard 
  4.         replace h_df`df' = h_df`df' * 1000
  5.         predict s_df`df', survival 
  6.         estimates store df`df'
=======
. 
. foreach df in 2 4 6 {
  2.         stpm2, df(`df') scale(hazard) bhazard(rate)
  3.         predict h_df`df', hazard 
  4.         replace h_df`df' = h_df`df' * 1000
  5.         predict s_df`df', survival 
  6.         estimates store df`df'
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
  7. }

Iteration 0:   log likelihood = -8857.8555  
Iteration 1:   log likelihood = -8602.0397  
Iteration 2:   log likelihood = -8598.8872  
Iteration 3:   log likelihood = -8598.8825  
Iteration 4:   log likelihood = -8598.8825  

Log likelihood = -8598.8825                              Number of obs = 7,775

------------------------------------------------------------------------------
             | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
xb           |
       _rcs1 |   .8587922    .027247    31.52   0.000      .805389    .9121954
       _rcs2 |   .2786853   .0239657    11.63   0.000     .2317133    .3256572
       _cons |  -1.820875   .0317465   -57.36   0.000    -1.883097   -1.758653
------------------------------------------------------------------------------
(7,775 real changes made)

Iteration 0:   log likelihood = -8743.3857  
Iteration 1:   log likelihood = -8588.5984  
Iteration 2:   log likelihood = -8588.1174  
Iteration 3:   log likelihood = -8588.1166  

Log likelihood = -8588.1166                              Number of obs = 7,775

------------------------------------------------------------------------------
             | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
xb           |
       _rcs1 |   .8217794   .0245014    33.54   0.000     .7737575    .8698014
       _rcs2 |   .2007402   .0225589     8.90   0.000     .1565255    .2449548
       _rcs3 |   .0786269    .012199     6.45   0.000     .0547173    .1025365
       _rcs4 |   .0018648   .0067393     0.28   0.782    -.0113439    .0150736
       _cons |  -1.812389   .0313923   -57.73   0.000    -1.873916   -1.750861
------------------------------------------------------------------------------
(7,775 real changes made)

Iteration 0:   log likelihood =  -8734.395  
Iteration 1:   log likelihood =  -8587.549  
Iteration 2:   log likelihood = -8587.1424  
Iteration 3:   log likelihood = -8587.1412  
Iteration 4:   log likelihood = -8587.1412  

Log likelihood = -8587.1412                              Number of obs = 7,775

------------------------------------------------------------------------------
             | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
xb           |
       _rcs1 |   .8216408   .0245061    33.53   0.000     .7736098    .8696719
       _rcs2 |   .1979202     .02278     8.69   0.000     .1532723    .2425681
       _rcs3 |   .0841377   .0133232     6.32   0.000     .0580248    .1102506
       _rcs4 |    .019898   .0084671     2.35   0.019     .0033027    .0364933
       _rcs5 |   .0010064   .0053114     0.19   0.850    -.0094038    .0114166
       _rcs6 |   .0059215   .0037762     1.57   0.117    -.0014797    .0133227
       _cons |  -1.812651   .0314002   -57.73   0.000    -1.874194   -1.751107
------------------------------------------------------------------------------
(7,775 real changes made)

. 
. /*Plot the excess hazard functions*/
. twoway (line h_df2 h_df4 h_df6 _t, sort lcolor(red blue black))

. 
. /*Plot the survival functions*/
. twoway (line s_df2 s_df4 s_df6 _t, sort lcolor(red blue black)) 

. 
. 
. estimates stats df2 df4 df6, n(2773)

<<<<<<< HEAD
Akaike's information criterion and Bayesian information criterion
=======
Akaike's information criterion and Bayesian information criterion
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

-----------------------------------------------------------------------------
       Model |          N   ll(null)  ll(model)      df        AIC        BIC
-------------+---------------------------------------------------------------
         df2 |      2,773          .  -8598.882       3   17203.76   17221.55
         df4 |      2,773          .  -8588.117       5   17186.23   17215.87
         df6 |      2,773          .  -8587.141       7   17188.28   17229.78
-----------------------------------------------------------------------------

. 

<<<<<<< HEAD

<<<<<<< HEAD

5. Fit a proportional excess hazards model

=======

5. Fit a proportional excess hazards model

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

%%stata 
stpm2 agegrp2 agegrp3 agegrp4 female year8594, bhazard(rate) ///
   df(4) scale(hazard) eform
predict h2, hazard per(1000) ci
predict s2, survival ci

<<<<<<< HEAD

. stpm2 agegrp2 agegrp3 agegrp4 female year8594, bhazard(rate) ///
=======
. stpm2 agegrp2 agegrp3 agegrp4 female year8594, bhazard(rate) ///
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
>    df(4) scale(hazard) eform

Iteration 0:   log likelihood = -8592.7514  
Iteration 1:   log likelihood = -8485.9814  
Iteration 2:   log likelihood = -8483.7774  
Iteration 3:   log likelihood = -8483.7688  
Iteration 4:   log likelihood = -8483.7688  

Log likelihood = -8483.7688                              Number of obs = 7,775

------------------------------------------------------------------------------
             |     exp(b)   Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
xb           |
     agegrp2 |   1.285544   .0963537     3.35   0.001     1.109911     1.48897
     agegrp3 |   1.729918   .1311554     7.23   0.000     1.491045    2.007059
     agegrp4 |   2.613592    .262209     9.58   0.000     2.147044     3.18152
      female |   .5819093   .0335882    -9.38   0.000      .519665    .6516092
    year8594 |   .6785354   .0390329    -6.74   0.000     .6061874     .759518
       _rcs1 |   2.308665   .0542546    35.60   0.000     2.204739    2.417489
       _rcs2 |   1.217353   .0260828     9.18   0.000      1.16729    1.269562
       _rcs3 |   1.084414   .0126098     6.97   0.000     1.059979    1.109413
       _rcs4 |   1.004209   .0065616     0.64   0.520     .9914306    1.017152
       _cons |   .1948568   .0131564   -24.22   0.000     .1707041    .2224268
------------------------------------------------------------------------------
Note: Estimates are transformed only in the first equation.

. predict h2, hazard per(1000) ci

. predict s2, survival ci

. 

<<<<<<< HEAD

6. Plot predicted excess hazard rate

=======

6. Plot predicted excess hazard rate

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

%%stata
twoway	(line h2 _t if agegrp1 == 1 & female == 0 & year8594 == 1, sort) ///
			(line h2 _t if agegrp4 == 1 & female == 0 & year8594 == 1, sort) /// 
<<<<<<< HEAD
			,legend(label (1 "Youngest") label (2 "Oldest")) scheme(sj) ///
			xtitle("Excess hazard rate") ytitle("Time since diagnosis")
=======
			,legend(label (1 "Youngest") label (2 "Oldest")) scheme(sj) ///
			ytitle("Excess hazard rate") xtitle("Time since diagnosis")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD

. twoway  (line h2 _t if agegrp1 == 1 & female == 0 & year8594 == 1, sort) ///
>                         (line h2 _t if agegrp4 == 1 & female == 0 & year8594 
> == 1, sort) /// 
>                         ,legend(label (1 "Youngest") label (2 "Oldest")) sche
> me(sj) ///
>                         xtitle("Excess hazard rate") ytitle("Time since diagn
> osis")
=======
. twoway  (line h2 _t if agegrp1 == 1 & female == 0 & year8594 == 1, sort) ///
>                         (line h2 _t if agegrp4 == 1 & female == 0 & year8594 
> == 1, sort) /// 
>                         ,legend(label (1 "Youngest") label (2 "Oldest")) sche
> me(sj) ///
>                         ytitle("Excess hazard rate") xtitle("Time since diagn
> osis")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

. 
. 

<<<<<<< HEAD

<<<<<<< HEAD

7. Time-dependent effects for age group

=======

7. Time-dependent effects for age group

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

%%stata 
stpm2 agegrp2 agegrp3 agegrp4 female year8594, bhazard(rate) df(3) scale(hazard) ///
			tvc(agegrp2 agegrp3 agegrp4) dftvc(2)
predict h3, hazard per(1000) ci
predict s3, survival ci

twoway	(line h3 _t if agegrp1 == 1 & female == 0 & year8594 == 1, sort) ///
			(line h3 _t if agegrp4 == 1 & female == 0 & year8594 == 1, sort) ///
<<<<<<< HEAD
			, legend(label (1 "Youngest") label (2 "Oldest")) scheme(sj) ///
			xtitle("Excess hazard rate") ytitle("Time since diagnosis")
=======
			, legend(label (1 "Youngest") label (2 "Oldest")) scheme(sj) ///
			ytitle("Excess hazard rate") xtitle("Time since diagnosis")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD

. stpm2 agegrp2 agegrp3 agegrp4 female year8594, bhazard(rate) df(3) scale(haza
=======
. stpm2 agegrp2 agegrp3 agegrp4 female year8594, bhazard(rate) df(3) scale(haza
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
> rd) ///
>                         tvc(agegrp2 agegrp3 agegrp4) dftvc(2)

Iteration 0:   log likelihood = -8586.4362  
Iteration 1:   log likelihood =  -8482.059  
Iteration 2:   log likelihood = -8479.6676  
Iteration 3:   log likelihood = -8479.6437  
Iteration 4:   log likelihood = -8479.6437  

Log likelihood = -8479.6437                              Number of obs = 7,775

------------------------------------------------------------------------------
             | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
xb           |
     agegrp2 |   .2741607   .0798607     3.43   0.001     .1176365    .4306849
     agegrp3 |    .555553   .0812003     6.84   0.000     .3964034    .7147026
     agegrp4 |    .934842    .110683     8.45   0.000     .7179073    1.151777
      female |  -.5457334   .0579363    -9.42   0.000    -.6592864   -.4321804
    year8594 |  -.3873942   .0576354    -6.72   0.000    -.5003576   -.2744309
       _rcs1 |    .851634   .0459294    18.54   0.000      .761614     .941654
       _rcs2 |   .1365924   .0357271     3.82   0.000     .0665685    .2066162
       _rcs3 |   .0697446   .0112343     6.21   0.000     .0477257    .0917635
_rcs_ageg~21 |  -.0210178   .0626366    -0.34   0.737    -.1437832    .1017477
_rcs_ageg~22 |   .0706612   .0480804     1.47   0.142    -.0235747     .164897
_rcs_ageg~31 |  -.0256869   .0665868    -0.39   0.700    -.1561946    .1048208
_rcs_ageg~32 |   .1174402   .0534022     2.20   0.028     .0127739    .2221065
_rcs_ageg~41 |   -.037214   .0856722    -0.43   0.664    -.2051285    .1307005
_rcs_ageg~42 |   .1407585   .0726802     1.94   0.053    -.0016921    .2832092
       _cons |  -1.655974   .0703345   -23.54   0.000    -1.793827   -1.518121
------------------------------------------------------------------------------

. predict h3, hazard per(1000) ci

. predict s3, survival ci

. 
. twoway  (line h3 _t if agegrp1 == 1 & female == 0 & year8594 == 1, sort) ///
>                         (line h3 _t if agegrp4 == 1 & female == 0 & year8594 
> == 1, sort) ///
<<<<<<< HEAD
>                         , legend(label (1 "Youngest") label (2 "Oldest")) sch
> eme(sj) ///
>                         xtitle("Excess hazard rate") ytitle("Time since diagn
> osis")
=======
>                         , legend(label (1 "Youngest") label (2 "Oldest")) sch
> eme(sj) ///
>                         ytitle("Excess hazard rate") xtitle("Time since diagn
> osis")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

. 

<<<<<<< HEAD

<<<<<<< HEAD

8. Predicted (time-dependent) excess hazard ratios

=======

8. Predicted (time-dependent) excess hazard ratios

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

%%stata

predict hr2, hrnum(agegrp2 1) ci
predict hr3, hrnum(agegrp3 1) ci
predict hr4, hrnum(agegrp4 1) ci

twoway	(line hr2 _t if agegrp2 == 1 & female == 0 & year8594 == 1, sort) ///
		(line hr3 _t if agegrp3 == 1 & female == 0 & year8594 == 1, sort) ///
		(line hr4 _t if agegrp4 == 1 & female == 0 & year8594 == 1, sort) ///
<<<<<<< HEAD
		, legend(label (1 "Agegrp2") label (2 "Agegrp3") label(3 "Agegrp4")) ///
		xtitle("Time since diagnosis") ytitle("Excess mortality rate") scheme(sj)
=======
		, legend(label (1 "Agegrp2") label (2 "Agegrp3") label(3 "Agegrp4")) ///
		xtitle("Time since diagnosis") ytitle("Excess mortality rate") scheme(sj)
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD

. 
=======
. 
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
. predict hr2, hrnum(agegrp2 1) ci

. predict hr3, hrnum(agegrp3 1) ci

. predict hr4, hrnum(agegrp4 1) ci

. 
. twoway  (line hr2 _t if agegrp2 == 1 & female == 0 & year8594 == 1, sort) ///
>                 (line hr3 _t if agegrp3 == 1 & female == 0 & year8594 == 1, s
> ort) ///
>                 (line hr4 _t if agegrp4 == 1 & female == 0 & year8594 == 1, s
> ort) ///
<<<<<<< HEAD
>                 , legend(label (1 "Agegrp2") label (2 "Agegrp3") label(3 "Age
> grp4")) ///
>                 xtitle("Time since diagnosis") ytitle("Excess mortality rate"
=======
>                 , legend(label (1 "Agegrp2") label (2 "Agegrp3") label(3 "Age
> grp4")) ///
>                 xtitle("Time since diagnosis") ytitle("Excess mortality rate"
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
> ) scheme(sj)

. 

<<<<<<< HEAD

%%stata
twoway	(rarea hr4_lci hr4_uci _t, sort pstyle(ci)) ///
		(line hr4 _t, sort), yline(1) legend(off) ///
<<<<<<< HEAD
        xtitle("Years from Diagnosis") ///
        ytitle("Excess Mortality Rate Ratio")
=======
        xtitle("Years from Diagnosis") ///
        ytitle("Excess Mortality Rate Ratio")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD

. twoway  (rarea hr4_lci hr4_uci _t, sort pstyle(ci)) ///
>                 (line hr4 _t, sort), yline(1) legend(off) ///
>         xtitle("Years from Diagnosis") ///
>         ytitle("Excess Mortality Rate Ratio")
=======
. twoway  (rarea hr4_lci hr4_uci _t, sort pstyle(ci)) ///
>                 (line hr4 _t, sort), yline(1) legend(off) ///
>         xtitle("Years from Diagnosis") ///
>         ytitle("Excess Mortality Rate Ratio")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

. 

<<<<<<< HEAD

<<<<<<< HEAD

9. Difference in relative survival

=======

9. Difference in relative survival

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

%%stata
predict sdiff4, sdiff1(agegrp4 1 female 0 year8594 1) ///
                sdiff2(agegrp4 0 female 0 year8594 1) ci
twoway  (rarea sdiff4_lci sdiff4_uci _t, sort pstyle(ci)) ///
        (line sdiff4 _t, sort), yline(0) legend(off) ///
<<<<<<< HEAD
        xtitle("Years from Diagnosis") ///
        ytitle("Difference in Relative Survival")
=======
        xtitle("Years from Diagnosis") ///
        ytitle("Difference in Relative Survival")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD

. predict sdiff4, sdiff1(agegrp4 1 female 0 year8594 1) ///
=======
. predict sdiff4, sdiff1(agegrp4 1 female 0 year8594 1) ///
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
>                 sdiff2(agegrp4 0 female 0 year8594 1) ci

. twoway  (rarea sdiff4_lci sdiff4_uci _t, sort pstyle(ci)) ///
>         (line sdiff4 _t, sort), yline(0) legend(off) ///
<<<<<<< HEAD
>         xtitle("Years from Diagnosis") ///
>         ytitle("Difference in Relative Survival")
=======
>         xtitle("Years from Diagnosis") ///
>         ytitle("Difference in Relative Survival")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

. 

<<<<<<< HEAD

<<<<<<< HEAD

10. Difference in excess mortality rate

=======

10. Difference in excess mortality rate

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

%%stata
predict hdiff4, hdiff1(agegrp4 1 female 0 year8594 1) ///
                hdiff2(agegrp4 0 female 0 year8594 1) ci
replace hdiff4 = hdiff4*1000
replace hdiff4_lci = hdiff4_lci*1000
replace hdiff4_uci = hdiff4_uci*1000
                
twoway  (rarea hdiff4_lci hdiff4_uci _t, sort pstyle(ci)) ///
        (line hdiff4 _t, sort), yline(0) legend(off) ///
<<<<<<< HEAD
        xtitle("Years from Diagnosis") ///
        ytitle("Difference in Excess Mortality Rate")
=======
        xtitle("Years from Diagnosis") ///
        ytitle("Difference in Excess Mortality Rate")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD

. predict hdiff4, hdiff1(agegrp4 1 female 0 year8594 1) ///
=======
. predict hdiff4, hdiff1(agegrp4 1 female 0 year8594 1) ///
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85
>                 hdiff2(agegrp4 0 female 0 year8594 1) ci

. replace hdiff4 = hdiff4*1000
(7,775 real changes made)

. replace hdiff4_lci = hdiff4_lci*1000
(7,775 real changes made)

. replace hdiff4_uci = hdiff4_uci*1000
(7,775 real changes made)

.                 
. twoway  (rarea hdiff4_lci hdiff4_uci _t, sort pstyle(ci)) ///
>         (line hdiff4 _t, sort), yline(0) legend(off) ///
<<<<<<< HEAD
>         xtitle("Years from Diagnosis") ///
>         ytitle("Difference in Excess Mortality Rate")
=======
>         xtitle("Years from Diagnosis") ///
>         ytitle("Difference in Excess Mortality Rate")
>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

. 

<<<<<<< HEAD

<<<<<<< HEAD

Muchas gracias por tu atencion. Thank you

=======

Muchas gracias por tu atencion. Thank you

>>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85

<<<<<<< HEAD ======= >>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85 <<<<<<< HEAD ======= >>>>>>> 38e07db0a280e7979c99143dc27e473f8ff73b85